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Abstract 

Most global optimization problems are nonlinear and thus difficult to solve, and they 
become even more challenging when uncertainties are present in objective functions and 
constraints. This paper provides a new two-stage hybrid search method, called Eagle 
Strategy, for stochastic optimization. This strategy intends to combine the random 
search using Levy walk with the firefly algorithm in an iterative manner. Numerical 
studies and results suggest that the proposed Eagle Strategy is very efficient for stochas- 
tic optimization. Finally practical implications and potential topics for further research 
will be discussed. 
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1 Introduction 



To find the solutions to any optimization problems, we can use either conventional 
optimization algorithms such as the Hill-climbing and simplex method, or heuristic 
methods such as genetic algorithms, or their proper combinations. Modern meta- 
heuristic algorithms are becoming powerful in solving global optimization problems 
[4, 6, 7, 9, 20, 21, 22], especially for the NP-hard problems such as the traveUing 
salesman problem. For example, particle swarm optimization (PSO) was developed 
by Kennedy and Eberhart in 1995 [8, 9], based on the swarm behaviour such as fish and 
bird schooling in nature. It has now been applied to find solutions for many optimization 
applications. Another example is the Firefly Algorithm developed by the first author 
[20, 21] which has demonstrated promising superiority over many other algorithms. The 
search strategies in these multi-agent algorithms are controlled randomization and ex- 
ploitation of the best solutions. However, such randomization typically uses a uniform 
distribution or Gaussian distribution. In fact, since the development of PSO, quite a 
few algorithms have been developed and they can outperform PSO in different ways 
[21, 23]. 
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On the other hand, there is always some uncertainty and noise associated with all 
real- world optimization problems. Subsequently, objective functions may have noise and 
constraints may also have random noise. In this standard optimization problem 

becomes a stochastic optimization problem. Methods that work well for standard opti- 
mization problems cannot directly be applied to stochastic optimization; otherwise, the 
obtained results are incorrect or even meaningless. Either the optimization problems 
have to be reformulated properly or the optimization algorithms should be modified 
accordingly, though in most cases we have to do both [3, 10, 19] 

In this paper, we intend to formulate a new metaheuristic search method, called 
Eagle Stategy (ES), which combines the Levy walk search with the Firefly Algorithm 
(FA). We will provide the comparison study of the ES with PSO and other relevant 
algorithms. We will first outline the basic ideas of the Eagle Strategy, then outline 
the essence of the firefiy algorithm, and finally carry out the comparison about the 
performance of these algorithms. 



2 Stochastic Multiobjective Optimization 

An ordinary optimization problem, without noise or uncertainty, can be written as 

min /i(x), (i = l,2,...,iV) (1) 

subject to (f>j{x) = 0, (j = 1, 2, J), 

V'fc(x) <0, (fc = l,2,...,i^), (2) 

where x — (xi,X2, ...,Xd)'^ is the vector of design variables. 

For stochastic optimization problems, the effect of uncertainty or noise on the design 
variable x, can be described by a random variable with a distribution Qi [10, 19]. 
That is 

Xi I-)- ^i{xi), (3) 

and 

6 - Qi- (4) 

The most widely used distribution is the Gaussian or normal distribution N{xi, Ui) with 
a mean Xi and a known standard deviation cji. Consequently, the objective functions 
fi become random variables /^(x, ^). 

Now we have to reformulate the optimization problem as the minimization of the 
mean of the objective function fi (x) or /z /. 

min{M/i,-,/^/„}. (5) 

Here = E{fi) is the moan or expectation of /,;(^(x)) where i = 1,2,...,N. More 
generally, we can also include their uncertainties, which leads to the minimization of 

min{/i/, + Xai, ...,/U/jv + Aajv}, (6) 

where A > is a constant. In addition, the constraints with uncertainty should be 
modified accordingly. 

In order to estimate Hf., we have to use some sampling techniques such as the Monte 
Carlo method. Once we have randomly drawn the samples, we have 



' p=l 



where Ni is the number of samples. 



3 Eagle Strategy 



The foraging behaviour of eagles such as golden eagles or Aquila Chrysaetos is inspiring. 
An eagle forages in its own territory by flying freely in a random manner much like the 
Levy flights. Once the prey is sighted, the eagle will change its search strategy to an 
intensive chasing tactics so as to catch the prey as efficiently as possible. There are two 
important components to an eagkrs hunting strategy: random search by Levy flight (or 
walk) and intensive chase by locking its aim on the target. 

Furthermore, various studies have shown that flight behaviour of many animals and 
insects has demonstrated the typical characteristics of Levy flights [5, 14, 12, 13]. A 
recent study by Reynolds and Frye shows that fruit flies or Drosophila melanogaster, 
explore their landscape using a series of straight flight paths punctuated by a sudden 
90° turn, leading to a Levy-flight-style intermittent scale-free search pattern. Studies 
on human behaviour such as the Ju/'hoansi hunter-gatherer foraging patterns also show 
the typical feature of Levy flights. Even light can be related to Levy flights [2]. Sub- 
sequently, such behaviour has been applied to optimization and optimal search, and 
preliminary results show its promising capability [12, 14, 16, 17]. 

3.1 Eagle Strategy 

Now let us idealize the two-stage strategy of an eagle's foraging behaviour. Firstly, we 
assume that an eagle will perform the Levy walk in the whole domain. Once it flnds a 
prey it changes to a chase strategy. Secondly, the chase strategy can be considered as 
an intensive local search using any optimization technique such as the steepest descent 
method, or the downhill simplex or Nelder-Mead method [11]. Obviously, we can also 
use any efficient metaheuristic algorithms such as the particle swarm optimization (PSO) 
and the Firefly Algorithm (FA) to do concentrated local search. The pseudo code of 
the proposed eagle strategy is outlined in Fig. 1. 

The size of the hypersphere depends on the landscape of the objective functions. If 
the objective functions are unimodal, then the size of the hypersphere can be about the 
same size of the domain. The global optimum can in principle be found from any initial 
guess. If the objective are multimodal, then the size of the hypersphere should be the 
typical size of the local modes. In reality, we do not know much about the landscape 
of the objective functions before we do the optimization, and we can cither start from 
a larger domain and shrink it down or use a smaller size and then gradually expand it. 

On the surface, the eagle strategy has some similarity with the random-restart hill 
climbing method, but there are two important differences. Firstly, ES is a two-stage 
strategy rather than a simple iterative method, and thus ES intends to combine a good 
randomization (diversiflcation) technique of global search with an intensive and efScient 
local search method. Secondly, ES uses Levy walk rather than simple randomization, 
which means that the global search space can be explored more efflciently. In fact, 
studies show that Levy walk is far more efficient than simple random-walk exploration. 

The Levy walk has a random step length being drawn from a Levy distribution 

Levy ~ u = t'^, (1< A < 3), (8) 

which has an infinite variance with an infinite mean. Here the steps of the eagle motion 
is essentially a random walk process with a power-law step-length distribution with a 
heavy tail. The special case A = 3 corresponds to Brownian motion, while A = 1 has 
a characteristics of stochastic tunneling, which may be more efficient in avoiding being 
trapped in local optima. 

For the local search, we can use any efficient optimization algorithm such as the 
downhill simplex (Nelder-Mead) or metaheuristic algorithms such as PSO and the firefiy 
algorithm. In this paper, we used the firefly algorithm to do the; local search, since the 
firefly algorithm was designed to solve multimodal global optimization problems [21], 



Eagle Strategy 

Objective functions /i(x), /jv(x) 
Initial guess x*"*^ 
while (1|x*+^ ~ x*|| > tolerance) 
Random search by performing Levy walk 
Evaluate the objective functions 
Intensive local search with a hypersphere 

via Nelder-Mead or the Firefly Algorithm 
if (a better solution is found) 

Update the current best 
end if 

Update t = t + 1 

Calculate means and standard deviations 
end while 

Postprocess results and visualization 

Figure 1: Pseudo code of the Eagle Strategy (ES). 



3.2 Firefly Algorithm 

We now briefly outline the main components of the Firefly Algorithm developed by 
the first author [20] , inspired by the flash pattern and characteristics of fireflies. For 
simplicity in describing the algorithm, we now use the following three idealized rules: 1) 
all fireflies arc unisex so that one firefly will be attracted to other fireflies regardless of 
their sex; 2) Attractiveness is proportional to their brightness, thus for any two flashing 
fireflies, the less brighter one will move towards the brighter one. The attractiveness 
is proportional to the brightness and they both decrease as their distance increases. If 
there is no brighter one than a particular firefiy, it will move randomly; 3) The brightness 
of a firefly is affected or determined by the landscape of the objective function. For a 
maximization problem, the brightness can simply be proportional to the value of the 
objective functions. 

In the firefly algorithm, there are two important issues: the variation of light inten- 
sity and formulation of the attractiveness. For simplicity, we can always assume that 
the attractiveness of a flrefly is determined by its brightness which in turn is associated 
with the encoded objective function. 

In the simplest case for maximum optimization problems, the brightnciss / of a firefiy 
at a particular location x can be chosen as J(x) oc /(x). However, the attractiveness /3 
is relative, it should be seen in the eyes of the beholder or judged by the other fireflies. 
Thus, it will vary with the distance rij between firefiy i and firefly j. In addition, light 
intensity decreases with the distance from its source, and light is also absorbed in the 
media, so we should allow the attractiveness to vary with the degree of absorption. In 
the simplest form, the light intensity I{r) varies according to the inverse square law 
J(r) = j± where Ig is the intensity at the source. For a given medium with a fixed light 
absorption coefficient 7, the light intensity I varies with the distance r. That is 

/ = /oe-^^ (9) 

where Iq is the original light intensity. 

As a firefly's attractiveness is proportional to the light intensity seen by adjacent 
firefiies, we can now define the attractiveness /3 of a firefly by 

/3 = /3oe-^'-', (10) 
where /3o is the attractiveness at r = 0. 



Firefly Algorithm 

Objective function /p(x), x = (xi, ...,Xd)'^ 
Initial population of fireflies Xj (i = 1, n) 
Light intensity li at is determined by /p(xi) 
Define light absorption coefficient 7 
while (t <MaxGeneration) 
for i = 1 : n all n fireflies 

for j = I : i all n fireflies 

if ft- > li) 

Move firefly i towards j (d- dimension) 

end if 

Vary /3 via exp[— 7r] 
Evaluate new solutions and update 
end for j 
end for i 

Rank the fireflies and find the current best 
end while 

Postprocess results and visualization 

Figure 2: Pseudo code of the firefly algorithm (FA). 



The distance between any two fireflies i and j at Xj and xj, respectively, is the 
Cartesian distance 



^ '^{Xi,k - Xj,k)^, (11) 
\ k=l 



where Xi^k is the A;th component of the spatial coordinate Xj of ith firefly. In the 2-D 
case, we have 

rij = ^Jixi-XjY + {yi-yjY. (12) 

The movement of a firefly i is attracted to another more attractive (brighter) firefly 
j is determined by 

2 1 
Xj = Xj + /3oe '^'''^ (xj - Xj) + a (rand - 2 )' (l'^) 

where the second term is due to the attraction. The third term is randomization with 
a control parameter a, which makes the exploration of the search space more efficient. 

We have tried to use different values of the parameters a, /3o,7 [20, 21], after some 
simulations, we concluded that we can use /3o = 1. a G [0,1], 7 = 1, and A = 1 for 
most applications. In addition, if the scales vary significantly in different dimensions 
such as —10^ to 10^ in one dimension while, say, —0.001 to 0.01 along the other, it is 
a good idea to replace a by aSk where the scaling parameters Sk(k = l,...,d) in the d 
dimensions should be determined by the actual scales of the problem of interest. 

There are two important limiting cases when 7 — > and 00. For 7^-0, the 
attractiveness is constant (3 = Po and the length scale T = 1/^7 00, this is equivalent 
to say that the light intensity does not decrease in an idealized sky. Thus, a flashing 
firefly can be seen anywhere in the domain. Thus, a single (usually global) optimum 
can easily be reached. This corresponds to a special case of particle swarm optimization 
(PSO) discussed earlier. Subsequently, the efficiency of this special case could be about 
the same as that of PSO. 

On the other hand, the limiting case 7 — >■ 00 leads to F — )• and /3(r) — ^ S{r) (the 
Dirac delta function), which means that the attractiveness is almost zero in the sight of 
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Figure 3: Ackley's function for two independent variables with a global minimum /* = at 



other fireflies or the fireflies are short-sighted. This is equivalent to the case where the 
fireflies roam in a very foggy region randomly. No other fireflies can be seen, and each 
firefly roams in a completely random way. Therefore, this corresponds to the completely 
random search method. As the firefly algorithm is usually in somewhere between these 
two extremes, it is possible to adjust the parameter 7 and a so that it can outperform 
both the random search and PSO. 

4 Simulations and Comparison 
4.1 Validation 

In order to validate the proposed algorithm, we have implemented it in Matlab. In our 
simulations, the values of the parameters are a = 0.2, 7=1,A = 1, and /3q = 1. As an 
example, we now use the ES to find the global optimum of the Ackley function 



where {d = 1, 2, ...) [1]. The global minimum /, = occurs at (0, 0, 0) in the domain 
of —32.768 < Xi < 32.768 where i = 1, 2, d. The landscape of the 2D Ackley function 
is shown in Fig. 3, while the landscape of this function with 2.5% noise is shown in Fig. 
4 

The global minimum in 2D for a given noise level of 2.5% can be found after about 
300 function evaluations (for 20 fireflies after 15 iterations, see Fig. 5). 



(0,0). 



/(x) = -20cxp[--. exp[-^cos(27ra;j)] + 20 + e, 




(14) 



4.2 Comparison of ES with PSO 



Various studies show that PSO algorithms can outperform genetic algorithms (GA) 
[7] and other conventional algorithms for solving many optimization problems. This is 



Figure 4: Ackley's 2D function with Gaussian noise. 



partially due to that fact that the broadcasting ability of the current best estimates 
gives better and quicker convergence towards the optimality. A general framework 
for evaluating statistical performance of evolutionary algorithms has been discussed in 
detail by Shilane et al. [15]. 

Now we will compare the Eagle Strategy with PSO for various standard test func- 
tions. After implementing these algorithms using Matlab, we have carried out extensive 
simulations and each algorithm has been run at least 100 times so as to carry out mean- 
ingful statistical analysis. The algorithms stop when the variations of function values 
are less than a given tolerance e < 10~^. The results are summarized in the following 
table (see Table 1) where the global optima are reached. The numbers are in the format: 
average number of evaluations (success rate), so 12. 7± 1.15(100) means that the average 
number (mean) of function evaluations is 12.7 x 10^ = 12700 with a standard deviation 
of 1.15 X 10^ — 1150. The success rate of finding the global optima for this algorithm is 
100%. Here we have used the following abbreviations: MWZ for Michalewicz's function 
with d = 16, RBK for Rosenbrock with d = 16, De Jong for De Jong's sphere function 
with d — 256, Schwefel for Schwefel with d = 128, Ackley for Ackley's function with 
d — 128, and Shubert for Shubert's function with 18 minima. In addition, all these test 
functions have a 2.5% of Gaussian noise, or ct = 0.025. In addition, we have used the 
population size n = 20 in all our simulations. 

We can see that the ES is noticeably more efficient in finding the global optima 
with the success rates of 100%. Each function evaluation is virtually instantaneous on 
a modern personal computer. For example, the computing time for 10,000 evaluations 
on a SGHz desktop is about 5 seconds. Even with graphics for displaying the locations 
of the particles and firefiies, it usually takes less than a few minutes. Furthermore, we 
have used various values of the population size n or the number of fireflies. We found 
that for most problems n = 15 to 50 would be sufficient. For tougher problems, larger 
n such as n = 100 or 250 can be used, though excessively large n should not be used 
unless there is no better alternative, as it is more computationally extensive. 
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Figure 5: The initial locations of the 20 fireflies (left) and their locations after 15 iterations 
(right). We have used 7 = 1. 



Table 1: Comparison of algorithm performance 





PSO (xlO^) 


ES (xlO^) 


Easom 


185.9 ±3.1(97) 


12.7 ± 1.15(100) 


MWZ 


346.1 ±8.0(98) 


36.1 ±3.5(100) 


Rosenbrock 


1637 ± 79(98) 


75 ± 6.4(100) 


De Jong 


852 ± 16(100) 


70.7 ± 7.3(100) 


Schwefel 


726.1 ±25(97) 


99 ±6.7(100) 


Ackley 


1170 ± 19(92) 


54 ±5.2(100) 


Rastrigin 


3973 ± 64(90) 


151 ± 14(100) 


Easom 


863.7 ± 55(90) 


76 ± 11(100) 


Griewank 


2798 ± 63(92) 


134 ±9.1(100) 


Shubert 


1197 ±56(92) 


32 ± 2.5(100) 



5 Conclusions 

By combining Levy walk with the firefly algorithm, we have successfully formulated a 
hybrid optimization algorithm, called Eagle Strategy, for stochastic optimization. After 
briefly outlining the basic procedure and its similarities and differences with particle 
swarm optimization, we then implemented and compared these algorithms. Our sim- 
ulation results for finding the global optima of various test functions suggest that ES 
can significantly outperform the PSO in terms of both efficiency and success rate. This 
implies that ES is potentially more powerful in solving NP-hard problems. 

However, we have not carried out sensitivity studies of the algorithm-dependent 
parameters such as the exponent A in Levy distribution and the light absorption coef- 
ficient 7, which may be fine-tuned to a specific problem. This can form an important 
research topic for further research. Furthermore, other local search algorithms such 
as the Newton-Raphson method, sequential quadratic programming and Nelder-Mead 
algorithms can be used to replace the firefly algorithm, and a comparison study should 
be carried out to evaluate their performance. It may also show interesting results if 
the level of uncertainty varies and it can be expected that the higher level of noise will 
make it more difficult to reach optimal solutions. 

As other important further studies, we can also focus on the applications of this 



hybrid algorithm on the NP-hard travehng salesman problem. In addition, many engi- 
neering design problems typically have to deal with intrinsic inhomogcneous materials 
properties and such uncertainty may often affect the design choice in practice. The 
application of the proposed hybrid algorithm in engineering design optimization may 
prove fruitful. 
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